###############################################
# Dynamic DiD: All Visas and Crime Rates
###############################################

rm(list=ls())

library(Hmisc)
library(readstata13)
library(foreign)
library(tidyverse)
library(did)
library(ggplot2)
set.seed(1111)

#sink("06_results_rates_all.txt")

############################
# Read and prepare data 
############################

# load data
load("final_county_2023april6.RData")
names(d_county)

#######################
# Results 
#######################

# d-did
out5 <- att_gt(yname = "crime_rates_s",
               gname = "first.treat_all",
               idname = "county",
               tname = "survey",
               data = d_county,
               est_method = "reg"
)

es5 <- aggte(out5, type = "dynamic",na.rm = TRUE)
summary(es5)
ggdid(es5)

# save results for overall effect
pe5 <- es5$overall.att
se5 <- es5$overall.se
results5 <- rbind(pe5,se5) 
results5
write.table(results5, "overallATT_5.txt", sep="\t")

# table formal to prepare plot
event_study_diff_never <-   data.frame(
  type          = "dynamic",
  term = paste0('ATT(', es5$egt, ")"),
  event.time= es5$egt,
  estimate  = es5$att.egt,
  std.error = es5$se.egt,
  conf.low  = es5$att.egt - es5$crit.val.egt * es5$se.egt,
  conf.high = es5$att.egt + es5$crit.val.egt  * es5$se.egt)

# plot
p_es_never3 <- ggplot(data = event_study_diff_never,
                      mapping = aes(x = event.time, y = estimate)) +
  geom_line(size = 0.5, alpha = 2, colour = "black") +
  theme_light()  +
  geom_vline(xintercept = 0-0.05, color = 'grey', size = 1.2, linetype = "dotted") + 
  geom_hline(yintercept = 0, colour="black",  linetype = "dotted")+
  geom_pointrange(aes(ymin = conf.low, ymax = conf.high), show.legend = FALSE, linetype= 1, size = 1.1,
                  color=ifelse(event_study_diff_never$event.time<0,"gray60","black"))+
  xlab("Survey wave since immigration shock") +
  ylab("DiD Estimate (in standard deviation units)") +
  scale_x_continuous(breaks = c(-5:3)) +
  scale_y_continuous(limits = c(-4,4)) +
  theme(axis.text.y = element_text(size = 14),axis.title.y = element_text(size = 18)) +
  theme(axis.text.x = element_text(size = 14),axis.title.x = element_text(size = 18)) +
  theme(axis.title = element_text(color="black",  size = 12))+
  theme(plot.title=ggtext::element_markdown(size = 12,
                                            color="black",
                                            hjust=0,
                                            lineheight=1.2)
  ) +
  annotate(geom="text", x=-0.5, y=3.5, label="Average of post-treatment coefficients (in standard deviation units): 
                                            0.038 
                                           (0.068)",
           color="black",
           size = 7)

p_es_never3

# figure 2 all
tiff('figure2_all.tiff', units="in", width=10, height=7, res=200, compression = 'lzw')
p_es_never3
dev.off()

# table a15
sink(file = "tableA15.txt")
summary(es5)
sink(file = NULL)

sink()
